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I. Grant Purpose 


The purpose of this Grant (NAG5-120) was to support a Guest Investi- 
gator study involving computer simulations of VLF wave-particle inter- 
actions observed in the magnetosphere from the ISEE-1 spacecraft. 


II. Period of Performance 
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The period of performance under this contract extended from October 1, 
1980, to September 30, 1981. 
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III. Work Provided 


The purpose of the present study was to develop a self-consistent 
computer code for particle simulations which could treat both Landau 
and cyclotron interactions. The main part of a two-and-one-half di- 
mensional (2-1/2D) electromagnetic (EM) simulation code (referred to 
as the 1 EM2 code' hereafter) was programmed during the time spent by 
Dr, H. Matsumoto as a Guest Investigator at Stanford University. The 
code treats the two-dimensional field structure as well as the three- 
dimensional particle motion that are Involved in interactions between 
obliquely propagating VLF waves and energetic electrons in the magne- 
tosphere. 

Figure 1 shows a sequence of operations in one time-step in the 
EM2 code where At is an increment of the time-step and n is an integer 
showing the step number. Four solid circles in the upper panel indi- 
cate the initial values given, at the n^* 1 time-step. Step 1 calculates 
the number density p on all grids by a scheme of area sharing of super--* 
particles with a shape of rectangular charge. In Step 2, we compute an 
electrostatic component of the electric field by solving the Poisson 
equation, An FFT technique is used in this calculation. In Step 3, 
we advance the magnetic field F for a half time-step by a difference 
equation corresponding to 

fh-vxF (1) 
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Step 1 Form p 

Step 2 Correct the Poisson Field 
Step 3 Advance B 

Step 4 Advance v by Buneman-Boris Method 

Step 5 Advance r 

Step 8 Form J 

Step 7 Advance r 

Step 8 Advance B 

Step 9 Advance E 


Figure 1. Leap-frog scheme and sequence of operations 
in one time-step in the EM2 Code. 


Since the magnetic and electric fields are updated to this time, 
the equation of motion of particles can now be solved in Step 4. We 
use a Buneman-Boris scheme for tne updating of the particle velocity 
7. Then we advance particle position r for the next half time-step 
in Step 5. As the particle velocities and positions are known, we 
next compute a current density 7 on each grid point by distributing 
a quantity q7 of each superparticle to the nearest four grid points 
by the area sharing scheme. In Step 7, the particle locations are 
updated. The magnetic field is then updated for the latter half 
time-step in Step 8. The electric field is finally updated in Step 
9 by a difference equation corresponding to 

||= c 2 (VxB-y o 0) (2) 

The correction of the electric field in the Poisson solver in Step 
2 is necessary and mandatory to compensate for the short-term violation 
of the equation of continuity during the leap-frog scheme for the par- f 
tide position 7 and velocity V. All difference equations for both the 
Maxwell equations and equation of motion are written in a time- and 
space-centered difference scheme which is time-reversible and of higher 
accuracy than forward or backward difference schemes. Figure 2 shows 
a staggered structure of the spatial grid points with full and half 
integer points on which the field, current, and charge quantities are 
assigned as indicated. This assignment and staggered leap-frog time 
scheme shown in Figure 2 are needed to construct the centered difference 
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Figure 2. Staggered grid system used in the EM2 Code. 


scheme for all equations, both in space and in time. Figure 3 shows the 
main program of the EM2 code. 
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The remaining subroutines and necessary diagnostic routines of the 
E?12 code were developed in Japan. The complete EM2 code requires more 
than 5,000 program statements. A detailed description of the code will 
be contained in a book entitled Computer Simulations in Space "Physios 
(edited by H. Matsumoto and T. Sato), which will be published from Terra- 
Pub/D. Reidel Book Publishing Company this autumn or next winter. 

In this report, we will present the results of two test runs of the 
EM2 code which demonstrate the validity of the code. One test run is 
designed to reproduce the dispersion characteristics of the normal modes 
of plasma slaves in a homogeneous magnetized plasma. The other is de- 
signed to reproduce the propagation characteristics of a monochromatic 
wave of a known mode. Figure 4 shows the result of the first test and 
Figure 5 gives those of the second test. Figure 4-1 shows the disper- 
sion of the measured wave fields after the test run, assuming an initial 
bi-Maxwell ian velocity with a temperature anisotropy of T a /T„ = 2. The* 
plasma frequency n e is set equal to 2 and cyclotron frequency is set 
equal to 1. All modes predicted by the linear theory of plasma waves 
are seen in the results. Figure 4-2 illustrates the configuration of 
wave magnetic and electric- field vectors as v/ell as of the wave number 
vector k for both 0- and X-mode waves in a homogeneous magnetized plasma. 
Figure 4-3 and 4-4 show the results of the dispersion characteristics 
(analyzed by an FFT technique) of the 0- and X-mode waves, respectively. 
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MAIN PROGRAM FOR WHISTLER INTERACTION BY 2-D EM CODE 
PROGRAMMED BY H. MATSUMOTO 
• JULY 1980 


1 COMMON/ 1 NPCQM/ 1 01 AG ,NSTEP .NQBNUM , ION 

2 COMMON/TIMECM/T 

3 COMMON/CONSTC/P, C, OX, DY, DT, NX, NY, NP, QS, DKX, DKY, AMI, 

NX1, NY1 , NPC, DXI , DYI, SLX, SLY, POT, ANPCI, 
P2DT, PDT2, PDTDX2, PDTDY2, PC2DTX, PC2DTY, 
QMEDT, QMIDT, NX2, NY2 



DATA INPUT 
START 
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INITIALIZATION 

CALL INPUT 

S 

MAIN LOOP 

CALL INITIAL 

6 

TIME 

DO 1000 I-l.MSTEP 

7 

JUDGE 

T - ( I-1+NSTEP*J0BNUM)*DT 

8 


J ■ M0D(I ,IDIAG) - 1 

9 

DIAGNOSTICS 

IF(J.NE.O) GO TO 100 

10 

B(N) TO B(N*0.5) 

CALL OIAGNS 

11 


CALL BFIELD 


CORRECTION OF E 


12 

B(N-0.5) TO B(N) 

CALL ECRRCT 

13 

V(N-0.5) TO V(N+0.5) 

CALL BFIELD 

14 

R(N) TO R( N+0 .5) 

CALL VELCTY 


Figure 3a. 


Part I of Main Program for VLF Wave Interaction. 
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15 

J(N+0.5) 

CALL POSITN 

16 

RCN+0.5) TO R(N+1) 

CALL CURRNT 

17 

RH0(N+1) 

CALL POSITN 

18 

E(N) TO E(N+1) 

CALL CHARGE 

19 


CALL EFIELD 

20 

1000 

STORE ALL DATA FOR 
CONSECUTIVE RUN 

CONTINUE 

21 


CALL STORE 

22 


STOP 

23 


END 


Figure 3b. Part II of Main Program of VLF Wave Interaction 
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Dispersion Check by EM2 Code 
Parallel Propagation 

( k // B 0 ) 
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Figure 4-1. Results of dispersion characteristics tests: 
Dispersion characteristics of waves propagating parallel 
to the external magnetic field. 
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Figure 4-2. Results of dispersion characteristics tests: 
Magnetic and electric field vectors along with k vector 
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Dispersion Relation by EM 2 CODE 
0 MODE .. 
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o A 1 Results of dispersion characteristics tests 
Dilp^rslo^ characteristics 'of O-.ode waves. 
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Dispersion Relation by EM2 CODE 
X MODE 
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Figur*? 4-4. Results of dispersion characteristics tests: 
Dispersion characteristics of X-mode waves. 
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An interesting feature of Figure 4-3 is that the stimulation results 
show not only the ordinary electromagnetic mode, which is predicted by 
the linear cold plasma theory, but it also shows the existence of ordi- 
nary electrogmagnetic cyclotron harmomic waves at multiple frequencies 
of the cyclotron frequency (fl e = 1). This mode was predicted from the 
kinetic linear theory by Russian physicists (e.g., Akhiezer et al . 
[1975]). Figure 4-4 shows the results for the X-mode analysis. In the 
B z and Ey component, the fast and slow (in other words, Z-mode) extra- 
ordinary modes are clearly seen. In the B z component, one can see a 
faint mode representing extraordinary electromagnetic cyclotron waves 
at the first and third multiples of the cyclotron frequency. Ir, the 
E x component, for which two figures are shown with small and large k 
regions, the electron Bernstein mode waves are clearly seen. All of 
these dispersion characteristics agree well with linear theory. Thus, 
as seen in Figure 4, the EM2 code reproduces all possible waves modes 
which are predicted by linear cold and kinetic theories. 

The second check of the code was performed by propagation tests p 
of three well-known modes. In Figure 5-1, the configuration of k 
vectors of those tested modes relative to the coordinate system and 
to the external magnetic field are shown. These three modes are a 
parallel propagating whistler wave, an electron plasma wave, and an 
obliquely propagating wave. The results of propagation tests of these 
waves are shown in Figures 5-2 to 5-4. All results show propagation 
characteristics which agree well with those predicted by the linear 
theory. In each figure, aj-k diagrams analyzed by FFTs are given with 
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theoretical solid lines calculated by the linear theory. 

Thus the development of the EM2 particle code has been completed 
and is now ready to be utilized for self-consistent simulation studies 
of nonlinear resonant interactions of electrons (and ions, if needed) 
with waves with k vector of arbitrary direction with respect to the 
external magnetic field, including a quantitative study of the ISEE 
wave and particle data. 
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ELECTRON PLASMA WAVE 



WHISTLER WAVE 

PARALLEL PROPAGATION 



WHISTLER WAVE 

OBLIQUE PROPAGATION 



Figure 5-1. Results of test wave propagation: Configuration 

of k and B 0 vectors. 
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WHISTLER WAVE 
PARALLEL PROPAGATION TEST 

(By Component) 
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Figure 5-2. Results of test wave propagation: Propagation 

and <o-k diagram of parallel propagating whistler wave. 
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ELECTRON PLASMA WAVE 
PROPAGATION TEST 

(Ex Component) 





n st 
X 
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Figure 5-3. Results of test wave propagation: Propagation 

and w-k diagram of electron plasma wave. 
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WHJSTL5R WAVE 
OBLIQUE PROPAGATION TEST 


time 


(By Component) 



Figure 5-4. Results of test wave propagation: Propagation 

and u-k diagram of obliquely propagating whistler wave. 
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